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We describe a program for computing the abundances of light elements produced during Big 
Bang Nucleosynthesis which is publicly available at http://parthenope.na.infn.it/. Starting from 
nuclear statistical equilibrium conditions the program solves the set of coupled ordinary differential 
equations, follows the departure from chemical equilibrium of nuclear species, and determines their 
asymptotic abundances as function of several input cosmological parameters as the baryon density, 
the number of effective neutrino, the value of cosmological constant and the neutrino chemical 
potential. The program requires commercial NAG library routines. 



Program summary 



Title of program: PArthENoPE 

Program URL: http://parthenope.na.infn.it/ 



Program obtainable from: parthenope@na.infn.it 

Computers: PC-compatible running Fortran on Unix, MS Windows or Linux 

Operating systems under which the program has been tested: Windows 2000, Windows XP, Linux 

Programming language used: Fortran 77 

External routines/libraries used: NAG libraries 

No. of lines in distributed program, including input card and test data: 4969 
No. of bytes in distributed program, including input card and test data: 192 Kb 
Distribution format: tar.gz 

Nature of physical problem: Computation of yields of light elements synthesized in the primordial 
universe. 

Method of solution: BDF method for the integration of the ODE's, implemented in a NAG routine 
Typical running time: 90 sec with default parameters on a Dual Xeon Processor 2.4GHz with 2.GB 
RAM 

PACS numbers: 26.35,+c DSF 13/07, FERMILAB-PUB-07-079-A, SLAC-PUB-12488 

"A I 'alta fantasia qui manco possa; 
ma gia volgeva il mio disio e 'I velle, 
si come rota ch'igualmente e mossa, 
I'amor che move il sole e V altre stelle." 

Dante Alighieri, "Commedia" - Paradiso, Canto XXXIII, 142-145 



I. INTRODUCTION 



Big Bang Nucleosynthesis (BBN) is one of the fundamental pillars of the Cosmological Standard Model. In the 
very early Universe, when the temperature of the primordial plasma decreased from a few MeV down to ~ lOkeV, 
light nuclides as 2 H, 3 He, 4 He and, to a smaller extent, 7 Li were produced via a network of nuclear processes. The 
relative abundances of these nuclear "ashes" with respect to hydrogen can be determined via several observational 
techniques and in different astrophysical environments. In the standard cosmological scenario and in the framework 
of the electroweak Standard Model, the dynamics of this phase is controlled by only one free parameter, the baryon 
to photon number density, which can thus be fixed by fitting experimental observations. This parameter can be 
also independently measured with very high precision by Cosmic Microwave Background anisotropies [l[ and the 
agreement with BBN result is quite remarkable. For reviews see e.g. Q in Q or 0]. 
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More in general, the quality and quantity of new cosmological and astrophysical data available in the last decade 
has led to an overall consistent picture of the evolution of the Universe, usually referred to as the "concordance" 
model. This is based on standard physics plus a few phenomenological parameters, for which an underlying theory 
is however still missing. At present, one thus faces the intriguing possibility that one might test models which 
go beyond our present understanding of fundamental interactions, in a way which is complementary to traditional 
earth-based laboratory and accelerator approaches. This is illustrated e.g. by the search for new light degrees of 
freedom which might contribute to the total energy density in the Universe in addition to photons and neutrinos. 
To pursue this programme, it is crucial to achieve a high level of accuracy in theoretical predictions for cosmological 
observables, at least at the level of experimental uncertainties. In the case of BBN, many steps have been done in 
this direction by a careful analysis of several key aspects of the physics involved in the phenomenon. The accuracy of 
the weak reactions which enter the neutron/proton chemical equilibrium has been pushed well below the percent level 
[E @> 0, fl @ • Similarly, the neutrino decoupling has been carefully studied by several authors by explicitly solving the 
corresponding kinetic equations [13, EH, E2, 13 EH EH • These two issues are mainly affecting the prediction of 4 He 
mass fraction, which presently has a very small uncertainty, of the order of 0.1 %, due to the experimental uncertainty 
on neutron lifetime. Finally, much study has been devoted to the analysis of several nuclear reaction rates entering 
the BBN network, as well as the corresponding uncertainties. This task involves a careful study of the available data 
or predictions on each reaction, including an update in light of new relevant experimental measurements, the choice 
of a reasonable protocol to combine them in order to obtain a best estimate and an error and, finally, the calculation 
of the corresponding thermal averaged rates. This issue has been extensively discussed in [9(, whose results have been 
used in the program described in the present paper, and (l6l. [l7j. An important benchmark in this development has 
been represented by the compilation of the NACRE Collaboration database [l8[ . 

In view of all these recent developments, we believe that the scientific community interested in BBN, in itself or as 
a tool to constrain new physics beyond the Standard Model, might find useful a new public BBN code which updates 
the pioneering achievements of [H, [H SH 1 . For this reason we have publicly released a code we have developed 
and continuously updated over almost a decade, which we named PArthENoPE and can be obtained at the URL 
|http:/ /parthenope.na.infn.it/| The aim of the present paper is to give a general description of the program and how 
to use it. After briefly summarizing in Section |TT] the theoretical framework of BBN and all major improvements 
implemented in PArthENoPE , we discuss in Section [TTTI a few extensions of the minimal standard BBN scenario which 
are also included in the code. In Section [IV] the main structures of PArthENoPE are outlined, while a comparison with 
the public code of [2fJ is discussed in Section [V] Finally, in Section |VT1 we report our conclusions. Hereafter we use 
natural units where the reduced Planck constant, the speed of light and the Boltzmann constant are fixed to 1, i.e. 
h = C = kg = 1. 



II. THE THEORY OF BIG BANG NUCLEOSYNTHESIS 



A. The set of equations 



We consider N nuc species of nuclides, whose number densities, 
density of baryons ns , 



rij, are normalized with respect to the total number 



X; 



II B 



n, p, 



2 H, 



(1) 



The list of all nuclides which are typically included in BBN analyses and considered in PArthENoPE is reported in 
Table D 

In the (photon) temperature range of interest for BBN, lOMeV > T > 0.01 MeV, electrons and positrons are 
kept in thermodynamical equilibrium with photons by fast electromagnetic interactions and distributed according to 
a Fermi-Dirac distribution function f e ± , with chemical potential fi e , parameterized in the following by the function 
(j) e = He/T . The pressure and energy density of the electro mag netic plasma (e ± and 7) is calculated in PArthENoPE by 
including the effect of finite temperature QED corrections [1 31 ] - Furthermore, electromagnetic and nuclear scatterings 
keep the non-relativistic baryons in kinetic equilibrium, and their energy density ps and pressure p^ are given by 



PB 



M u + Y, ( AAf, + -T I X. 



n B , 



(2) 



1 It is a pleasure to acknowledge the public code of [2dl | as the starting point for many scholars interested in BBN, including the authors 
of the present paper. 
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No. 


Nuclide 


No. 


Nuclide 


No. 


Nuclide 


No. 


Nuclide 


No. 


Nuclide 


1 


n 


7 


6 Li 


13 


io B 


19 


13 C 


25 


15q 


2 


P 


8 


7 Li 


14 


n B 


20 


13 N 


26 


16 Q 


3 


2 H 


9 


7 Be 


15 


"C 


21 


14 C 






4 


3 H 


10 


8 Li 


16 


12g 


22 


14 N 






5 


3 He 


11 


8 B 


17 


12 C 


23 


14q 






6 


4 He 


12 


9 Be 


18 


12 N 


24 


15 N 







TABLE I: Nuclides considered in PArthENoPE 



Pfl 



= Tn B > (3) 



with AMj and M u the i-th nuclide mass excess and the atomic mass unit, respectively. 

The set of differential equations ruling primordial nucleosynthesis is the following (see for example @, BE^j): 



d 8ttG n 

~a =H= V— p ' (4) 

^ = -3H , (5) 
n B 



p=-3H(p + p) , (6) 

/ x n, x N k x m x N o \ 

x i = E N * y Tki ^ ~ Ti ^ kl wnIi )^ Ti ' (7) 

Y,Z 3 X 3 =n e --n e+ ^L(^^ e )^T"L{^^ e ) , (8) 



j,k,l 

n B 



where p and p denote the total energy density and pressure, respectively, 

P = Pj + Pe + Pv + Pb , (9) 
p = p 7 + p e + p„ + ps , (10) 

while i,j,k,l denote nuclear species, Ni the number of nuclides of type i entering a given reaction (and analogously 
Nj, Nk, Ni), and the T's denote symbolically the reaction rates. For example, in the case of decay of the species i, 
Ni — 1, Nj — and J2^i^kl is the inverse lifetime of the nucleus i; for binary collisions, N = Nj = Nk = N — 1 
and Tij^ki = {&ij^ki v), i.e. it represents the thermal average of the cross section for the reaction i + j ^ k + I times 
the relative velocity of i and j. In Eq. (8), Zi is the charge number of the i— th nuclide, and the function L(£,u>) is 
defined as 



i r°° 



i i 



(11) 



Equation ^ is the definition of the Hubble parameter H, a denoting the scale factor of the Friedmann-Robertson- 
Walker-Lemaitre metric, with Gat the gravitational constant, whereas Eq.s ([5]) and ([6]) state the total baryon number 
and entropy conservation per comoving volume, respectively. The set of N nuc Boltzmann equations ([7]) describes the 
density evolution of each nuclide specie, with Tki^ij the rate per incoming particles averaged over kinetic equilibrium 
distribution functions. Finally, Eq. © states the Universe charge neutrality in terms of the electron chemical potential, 
with L (m e /T, <p e ) the charge density in the lepton sector in unit of the electron charge. 

The neutrino energy density and pressure are defined in terms of their distributions in momentum space as 

Pv = 3p„ = 2 J ^ \p\ [f Ve + 2 /„.] , (12) 

Indeed, in the default scenario we assume a vanishing neutrino chemical potential, so that f Ue — fp e and f„ x = f v — 
fa = f Vr = f,j T . The nuclide evolution can be followed in PArthENoPE also for finite neutrino chemical potential, see 
Section fllT] below. 
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As well known, neutrinos decouple from the electromagnetic plasma at temperatures of a few MeV. Soon after, when 
the onset of e + — e~ annihilations takes place, are still partially coupled to neutrinos. The neutrino distributions 
are thus slightly distorted, especially in their high energy tail (and the e— flavor more than the other two, since the 
former also interacts via charged current). To get BBN predictions accurate at the sub-percent level it is necessary 
to follow in details this residual out of equilibrium neutrino heating by solving the kinetic equations for neutrino 
distributions. Remarkably, baryons provide a negligible contribution to the dynamics of the Universe at the BBN 
epoch as the baryon to photon number density is very small, rj < 10~ 9 , and therefore Boltzmann equations for neutrino 
species can be solved together with equations Q and ^ only, ignoring the dynamics of nuclear species. This allows 
one to solve the evolution of the neutrino species first, and then to substitute the resulting neutrino distribution into 
the remaining equations. We do not consider neutrino oscillations, whose effect has been studied and shown to be 
sub- leading in [Tj] ■ The reader can find further details on the neutrino decoupling stage in [§, [TJ, [l5[ . 



B. Numerical solution of the BBN set of equations 

The BBN set of equations (H])-© can be recast in a form more convenient for a numerical solution, which follows 
the evolution of the AT nuc + 1 unknown quantities (<p e , Xj) as functions of the dimensionless variable z — m e /T. In 
this framework, Eq. JS]) provides % as a function of <f) e . In particular, the set of differential equations implemented 
in PArthENoPE is the following: 



l Lki+ (p el B + Pe 7 S + J «2 

(pe-,B + Pe 7 S + ^T^) 



dz Z f Msl 9L 



(13) 



where 



dz l dz 3z h IMs.-ai L R+d R + msl) ' [ ' 



V Zifi 

dL - 

K 2 = z — 31/ — Z L— =r . (16) 



d 



z 



3 

3 



See Appendix fAI for notations and the explicit derivation of this set of equations. Equations (fl3|) and (fT4]) are solved 
by imposing the following initial conditions at z,„ = m e /(10 MeV): 

<t>e{Zin) = , (17) 

Xi(z in ) = X n (z m ) = (exp{qz in } + ly 1 , (18) 
X 2 (z in ) = X p (z in ) = (exp{-qz in } + 1) _1 , (19) 



(z m ) exp {b, zin} i = 2 H, 3 H, ... . (20) 

In the previous equations q — (M n —M p )/m e , and the quantities A, and Bi denote the atomic number and the binding 
energy of the i— th nuclide normalized to electron mass, respectively. Also note that Eq. (20) is only applied if the 
resulting abundance is greater than the numerical zero assumed (variable YMIN, whose default setting is 10~ 30 ). 
Finally, rji is the initial value of the baryon to photon number density ratio at T = 10 MeV (for a discussion of how it 
is related to the final value after e + — e~ annihilation stage see e.g. Section 4.2.2 in 0), and (f> e the solution of the 
implicit equation 

L(z m , <f> e °) - m £ Zi Xi(zin) ■ (21) 
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No. 


Reaction 


Type 


No. 


Reaction 


Type 


1 


m->p 


weak 


22 


6 Li 4- p «- 


■+ 7 + 7 Be 


(P,7) 


2 


3 H -> v e + e~ + 3 He 


weak 


23 


6 Li + p <-> 


3 He 4- 4 He 


3 He Pickup 


3 


8 Li v e + e- +2 4 He 


weak 


24 


7 Li 4- p <-> 


4 He 4- 4 He 


4 He Pickup 


4 


12 B -> ^ e + e- + 12 C 


weak 


24 bis 


7 Li 4- p <-> 7 


+ 4 He + 4 Hc 


(P7) 


5 


14 C -> P e + e" + 14 N 


weak 


25 


4 He 4- 2 H 


<-> 7 4- 6 Li 


(d,7) 


6 


S B -> z/ e + e+ + 2 4 He 


weak 


26 


4 He 4- 3 H 


^ 7 + 7 Li 


(t,7) 


7 


n C ^u e + e+ + "B 


weak 


27 


4 He 4- 3 He 


«74 7 Be 


( 3 He, 7 ) 


8 


12 N -> i/ e + e+ + 12 C 


weak 


28 


2 H 4- 2 H < 


-> n + 3 He 


2 H Strip. 


9 


13 N ^ u e + e + + 13 C 


weak 


29 


2 H + 2 H 


<-> p + 3 H 


2 H Strip. 


10 


14 - u e + e+ + 14 N 


weak 


30 


3 H 4- 2 H < 


-» n + 4 He 


2 H Strip. 


11 


ls O -> v e + e+ + 15 N 


weak 


31 


3 He 4- 2 H 


^ p + 4 He 


2 H Strip. 


12 


p 4- n «-» 7 + 2 H 


0,7) 


32 


3 He 4- 3 He <-i 


• p 4- p 4- 4 He 


( 3 He,2p) 


13 


2 H + n <-> 7 + 3 H 


0,7) 


33 


7 Li 4- 2 H <-> n 4- 4 He + 4 He 


(d,n a) 


14 


3 He + n ^ 7 + 4 He 


0,7) 


34 


7 Be 4- 2 H <-> p 4- 4 He + 4 He 


(d,p a) 


15 


6 Li 4- n <-> 7 4- 7 Li 


0*7) 


35 


3 He 4- 3 H 


«-» 7 4- 6 Li 


(t,7) 


16 


3 He + n« p + 3 H 


charge ex. 


36 


6 Li 4- 2 H • 


*-» n 4- 7 Be 


2 H Strip. 


17 


7 Be 4- n <-> p 4- 7 Li 


charge ex. 


37 


6 Li + 2 H 


<-> p 4- 7 Li 


2 H Strip. 


18 


6 Li 4- n <-> 3 H 4- 4 He 


3 H Pickup 


38 


3 He + 3 H ♦ 


-» 2 H 4- 4 He 


( 3 H,d) 


19 


7 Bc + ri <-> 4 He + 4 He 


4 He Pickup 


39 


3 H + 3 H <-> 


n 4- n 4- 4 He 


(t,nn) 


20 


2 H + p <-> 7 4- 3 He 


(P,7) 


40 


3 He + 3 H <-> 


p 4- n + 4 Hc 


(t,n p) 


21 


3 H4p«74 4 He 


(P»7) 











TABLE II: The reactions used in the small network. 
C. The Nuclear Chain 

In Tables [Til Hill and IIVI are reported the nuclear processes considered in PArthENoPE . The enumeration shown in 
the first column of the tables correspond to the order in which they appear in the program. See Q for the relevant 
formalism concerning the thermally averaged nuclear rates and an analysis of the main experimental reaction rates. 
Reactions included in Table IIT1 are used when running PArthENoPE in its simpler version (small network), while those 
of Tables IIIII and IIVI are added in the intermediate and complete network running options, which also follows the 
evolution of the nuclides heavier than 7 Be and 12 N, respectively. Using the small network gives values of the lighter 
nuclides like 2 H, 3 He, 4 He and 7 Li which differ from the results obtained with the complete network for less than 
0.02%, for default values of the input cosmological parameters. With respect to the database used in there 
are a few minor upgrades implemented here, namely the three reactions (98, 99, 100) have been inserted following 
the analysis of the extended network reported in |22j |. Also, we have added the recent data reported in [23| to the 
regressions for the rates (28, 29), with results in good agreement with those adopted in Q. 

III. NON-STANDARD PHYSICS 

In the standard scenario the only free parameter entering the BBN dynamics is the value of the baryon to photon 
number density 77, or equivalently the baryon energy density parameter fish 2 , see e.g. Q for the relation between 
these parameters. If one goes beyond the standard framework, the B BN p redictions may be altered by non-standard 
physics entering e.g. the neutrino [3, EH HH HE] or gravity sector |26|.|27||. o r more generically by the presence in the 
plasma of other degrees of freedom besides the Standard Model ones Hfli |U [H US, H SI . For an earlier review, 
see [35J. Typically, constraints on non-minimal and/or exotic scenarios require model-dependent modifications of the 
equations ruling BBN. In PArthENoPE we implement a few of them described below which are general enough to be 
commonly used/referred to in the specialized literature. 
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No. 


Reaction 


Type 


No. 


Reaction 


Type 


41 


7 Li + 3 H <- 


-» n + 9 Be 


3 H Strip. 


58 


6 Li + 4 He <-> 7 + 10 B 


(«,7) 


42 


7 Be + 3 H < 


-> p + 9 Be 


3 H Strip. 


59 


7 Li + 4 He <-> 7 + n B 


(a,7) 


43 


7 Li + 3 He - 


-» p + 9 Be 


3 He Strip. 


60 


7 Be + 4 He ^ 7 + n C 


(a,7) 


44 


7 Li + n <- 


+ 7 + 8 Li 


(n,7) 


61 


8 B + 4 He <-> p + U C 


(a.P) 


45 


10 B + n «- 


+ 7 + U B 


( n :7) 


62 


8 Li + 4 He <-> n + n B 


(a,n) 


46 


ll B + n«- 


* 7 + 12 B 


( n .7) 


63 


9 Be + 4 He <-> n + 12 C 


(a,n) 


47 


u C + n<- 


♦ p + u B 


(n,p) 


64 


9 Be + 2 H «-> n + 1() B 


( 2 H,n) 


48 


10 B + n «-» 


4 He + 7 Li 


(n,a) 


65 


10 B + 2 H <-> p + n B 


( 2 H,P) 


49 


7 Be + p <■ 


-» 7 + 8 B 


(PjT) 


66 


n B + 2 H«n+ 12 C 


( 2 H,n) 


•50 


9 Be + p ^- 


* 7 + 10 B 


(P.7) 


67 


4 He + 4 He + n <-> 7 + 9 Be 


(an,7) 


51 


10 B + p <- 


♦ 7 + U C 


(P»7) 


68 


4 He + 4 He + 4 He <-> 7 + 12 C 


(qq,7) 


52 


"B + p^ 


* 7 + 12 C 


(P>7) 


69 


8 Li + p <-> n + 4 He + 4 He 


(p,na) 


53 


u C + p«7 + 12 N 


(P.7) 


70 


8 B + n <-> p + 4 Hc + 4 Hc 


(n,pa) 


54 


12 B + p^ 


♦ n + 12 C 


(P. n ) 


71 


9 Be + p ^ 2 H + 4 He + 4 He 


(p,da) 


55 


9 Be + p <-> 


4 He + 6 Li 


(p.a) 


72 


U B + p «-> 4 He + 4 He + 4 He 


(p,aa) 


56 


10 B + p ^ 


4 He + 7 Be 


(p.a) 


73 


"C + n «■ 4 He + 4 He + 4 He 


(n,a a) 


57 


12 B + p ^ 


4 He + 9 Be 


(p. q ) 









TABLE III: The reactions used in the intermediate network in addition to those of Table [TT] 



A. Energy density of the vacuum, pa 

As in the original Kawano code we allow for a non-zero cosmological constant term at the BBN epoch. We 
parameterize it by means of p\ entering the equations only via 



3H -> BH 



'24 7tG 



N 



/m e \ 4 

{—) P + PA 



(22) 



The allowed range for this parameter in units of MeV 4 is < (pA/MeV 4 ) < 1. 



No. Reaction Type No. Reaction Type 



74 


12 C + n 


^ 7 + 13 c 


( n .7) 


88 


12 c + 


4 He 


<-» 7 + 16 


(«,7) 


75 


13 C + n 


«- 7 + 14 c 


( n .7) 


89 


10 B + 


4 He 


^ P + 13 c 


(a.P) 


76 


14 N + n 


<-> 7 + 15 N 


( n >7) 


90 


n B + 


4 He 


«-» P + 14 c 


(a.p) 


77 


13 N + n 


^ P + 13 C 


(n.P) 


91 


n C + 


4 He 


<-> p + 14 N 


(a.p) 


78 


14 N + n 


^ P + 14 c 


( n .P) 


92 


12 N + 


4 He 


^ P + 15 o 


(«>P) 


79 


ls O + n 


^ p + 15 N 


( n .P) 


93 


13 N + 


4 He 


<-> p + ie o 


(«>P) 


80 


15 + n «• 


-> 4 He + 12 C 


(n,a) 


94 


10 B + 


4 He 


«-» n + 13 N 


(a,n) 


81 


12 C + p 


<^ 7 + 13 N 


(p.7) 


95 


14 B + 


4 He 


<-> n + 14 N 


(a,n) 


82 


13 C + p 


<-> 7 + 14 N 


(p.7) 


96 


12 B + 


4 He 


»n + 15 N 


(a,n) 


83 


14 C + p 


<-» 7 + 15 N 


(p>7) 


97 


13 C + 


4 He 


«-> n + 16 


(a,n) 


84 


13 N + p 


«-> 7 + 14 


(P.7) 


98 


n B + 


2 H 


>-> P + 12 B 


2 H Strip. 


85 


14 N + p 


<-> 7 + 15 


(P.7) 


99 


12 C + 


2 H 


- P + 13 C 


2 H Strip. 


86 


15 N + p 


<-> 7 + 1G 


(P.7) 


100 


13 c + 


2 H 


^ P + 14 c 


2 H Strip. 


87 


15 N + p <- 


-> 4 He + 12 C 


(p,a) 













TABLE IV: The reactions used in the complete network in addition to those of Tables HTl and ITTTl 
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B. Extra degrees of freedom, AN c tc 



We parameterize the radiation density in non-elcctromagnetically interacting particles at the BBN epoch by an 
additional radiation energy density px entering H. This is related to the "number of extra effective neutrino species" 
customarily used in the literature AN C ^ by the equation 



Px 



8 30 



AiVeffTl , 



(23) 



where, from the entropy conservation, Tx = T = m e /z at temperatures higher than the effective neutrino decoupling 
temperature, chosen as Td — 2.3 MeV, or else 



T. 



x 



Pe,j,B(T) + p e , 7i s(T) 



"I 1/3 



_Pe,f,B(Td) + p e , 7 ,s(Td) 

The user may input a value of AiV ff in the range —3.0 < AN c g < 15.0 



, T <T d 



(24) 



C. Chemical potential of the neutrinos, £ 



The usual argument in favor of a cosmic lepton asymmetry is that sphaleron effects before electroweak symmetry 
breaking equilibrate the lepton and baryon asymmetries to within a factor of order unity, thus producing the ob- 
served baryon density. In principle, however, the electron-neutrino degeneracy parameter £ = [i Vs ,jT Va as well as the 
degeneracy parameters of the other neutrino flavors are not determined within the Standard Model, and should be 
constrained observationally. Recently, it has been realized that the measured neutrino oscillation parameters imply 
that neutrinos reach approximate chemical equilibrium before the BBN epoch. Thus, all neutrino chemical potentials 
can be taken to be equal, i.e. they are all characterized by the same degeneracy parameter £ that applies to v e 
[III H3j HI] • In light of these results it is meaningful to assume a single and shared value £ as the only free input 
parameter. Also, to achieve an approximate agreement between the observed and predicted light element abundances 
a possible lepton asymmetry must be small, |£| <C 1. For such small £ values the most important impact on BBN is 
a shift of the beta equilibrium between protons and neutrons. A subleading effect is a modification of the radiation 
density, 

^ 3 [f(!) 2 + H(iy : 

Moreover, the neutrino decoupling temperature is higher than in the standard case [39|, |40j, so that in principle one 
could get a non-standard T U (T) evolution, but such effects are completely negligible for our case. A non-zero £ slightly 
modifies the partial neutrino reheating following the e + e _ annihilation [4l| . again a completely negligible effect for 
the range of £ of our interest. In the code, we allow the user to select among 21 possible values of £, between -1.0 and 
+1.0, spaced by 0.1. The changes in the weak reactions are then automatically implemented, as in [8]. The associated 
change in AN e g of Eq. (j25jl is also accounted for. Note that to derive results on a finer grid a perturbative approach 
as the one in (3~H would be required. This possibility is left for an implementation in a future upgraded version of the 
code. 



(25) 



IV. THE STRUCTURE OF PARTHEN0PE 



The code is divided in two files, main.f and parthenope . f, the former one containing the main program and the 
latter the remaining subroutines. While all the physics is implemented in parthenope . f , the file main . f is an interface 
which can be possibly adapted to the user needs. The user can choose between two running modes: an interactive 
one, with parameter selections given on the screen, and a card mode requiring an input card, an example of which 
is provided as the file input (see also Table IVij) . The program links to the NAG libraries [42[ for some algebraic 
operations and the evaluation of special functions. 

The logical structure of the code is depicted in Figure [U In the following we detail the aim of each block. 
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•IfTEHD 




OUTEUOL 



FIG. 1: The logical structure of PArthENoPE 



MAIN 



MAIN contains the interface which allows the user to choose the physical and network input parameters and to 
customize the output. 

Physical parameters presently are: baryon density, number of additional neutrino species, neutron lifetime, neutrino 
chemical potential, energy density of the vacuum at the BBN epoch. 

Network parameters include: the choice among a small (9 nuclides and 40 reactions), an intermediate (18 nuclides 
and 73 reactions), and a complete network (26 nuclides and 100 reactions). Moreover, the user can change the rates 
of each reaction included in the chosen network, selecting a 'LOW or a 'HIGH' value, based on the experimental or 
theoretical uncertainties, or a customized multiplicative 'FACTOR'. 

Finally, the output options include: the choice of the nuclides whose evolution has to be followed versus z, the name 
of the output files (a first one with the final results and a second with the evolution of the selected nuclides), and the 
possibility to follow the status of the evolution on the screen. 



KEYWORD DESCRIPTION DEFAULT RANGE/OPTIONS 



OMEGABH 


Baryon density Qsh 2 


0.0223 


0.01 ~ 0.03 


DNNU 


Number of additional neutrino species 


0. 


-3. -r- 15. 


TAU 


Neutron lifetime 


885.7 s 


880. ~ 890. 


IXIE 


Integer fixing the electron neutrino chemical potential 


11 


1 -r- 21 


RHOLMBD 


Energy density of the vacuum, /9a, in MeV 4 


0. 


0. -f- 1. 


NETWORK 


Number of nuclides in the network 


9 


9,18,26 


FOLLOW 


Option for following the evolution on the screen 


F 


T,F 


OVERWRITE 


Option for overwriting the output files 


F 


T,F 


FILES 


Name of the output files 


parthenope.out 


20 bit string 






nuclides. out 


20 bit string 


OUTPUT 


Evolution of nuclides 


first 9 nuclides 


see text 


RATES 


Details on changed rates 


No change 


see text 



EXIT Closing keyword in the input card 



TABLE V: The list of the possible keywords in the input card, their default values and corresponding ranges/options. 
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RATES 


3 ( 12 1 0. ) ( 28 3 .4) (29 2 ) 


options for changing the nuclear rates 


RATES 


2 ( 3 2 0. ) (5 3 .6) 


options for changing the nuclear rates 


TAU 


885.7 


experimental value of neutron lifetime 


DNNU 


.0 


number of extra neutrinos 


IXIE 


11 


integer giving the value of v R chemical potential 


RHOLMBD 


.0 


value of cosmological constant energy density at the BBN epoch 


OVERWRITE 


T 


option for overwriting the output files 


FOLLOW 


T 


option for following the evolution on the screen 


OMEGABH 


.0223 


value of Q,sh 2 


NETWORK 


9 


number of nuclides in the network 


FILES 


parthenopel.out nuclidesl.out 


names of the two output files 


OUTPUT 


T 3 2 3 4 


options for customizing the output 


EXIT 




terminates input 




TABLE VI: 


An example of input card. 



All this information can be provided either interactively, following the on-screen instructions, or by an input card, 
with the format of the example card included in the distribution. In particular, each line in this card must start with 
an allowed key and the last line key has to be 'EXIT'. In order to be recognized, each key must start at the first bit of 
the line. The allowed keywords are listed in Table [V] together with the default values adopted by the code whenever 
the corresponding key is not explicitly set. 

An example of input card is shown in Table EO While the keys OMEGABH, DNNU, TAU, IXIE, RHOLMBD, 
NETWORK, FOLLOW, and OVERWRITE have only one argument, the keyword FILES has two arguments, that 
is the two names of the output files, each one at most 20 bits long. Some more details deserve the two keywords 
OUTPUT and RATES. OUTPUT can have at most N nuc + 2 arguments {N nuc being the number of nuclides of the 
chosen network), which are: 1) a bit equal to 'T' or 'F', if the user wants or not to store in the output the evolution 
of a given set of nuclide abundances, 2) the total number of such nuclides, and 3) the identity of these nuclides (given 
as the corresponding number in Table [J). In the example of Table IVTl with the sequence T 3 2 3 4 the user has chosen 
to store in the output the three nuclides p, 2 H and 3 H. The keyword RATES allows to change the default values of 
the nuclear rates used in the chosen network. The input card can have more than one line with this keyword, as in 
the example of Table IVTl Each line contains: 1) an integer k, giving the number of reactions whose change is specified 
on that line; 2) the kind of change for the k reactions with the syntax (m i f ), indicating that the reaction number 
m (see Tables HI1 IIIII and HVl) has to be changed according to the type of change i, with the factor f . Obviously, m can 
assume the values 1,...,M (M being the number of reactions of the chosen network), while i=l,2,3 corresponds to the 
'LOW', 'HIGH', or 'FACTOR' type of change, respectively. Whenever a statistically sound analysis is possible, as it 
is the case for most of the main reactions, the 'LOW'/'HIGH' rates represent 1 u lower/upper limits to the rate, as 
compiled in [9J. For most of the subleading reactions, they represent estimated ranges of variability, obtained from the 
literature. The main (or unique) reference from which the rate has been taken is reported as a comment in the code 
next to the related reaction line. Finally, if i=3 the real number f is the value of the multiplicative factor applied to 
the chosen reaction rate (not considered if the options i=l,2 are selected). For example, the first line of the input 
card of Table IVII specifies that 3 reaction rates should be changed in running PArthENoPE as follows 

p + n <-» 7 + 2 H low rate 
2 H + 2 H <-> n + 3 He rate multiplied by the factor 0.4 
2 H + 2 H <-> p + 3 H high rate 

Notice that it is possible to add comments after the parameters in the input card and the order of the lines with 
different keywords is not important. 

B. PARTHENOPE 

This subroutine drives the resolution of the BBN set of equations. It starts calling the initialization routine INIT, 
then the NAG solver, finally the output printing routine OUTEND. The NAG resolution parameters, controlling for 
example the resolution method and the numerical accuracy, have been chosen to optimize the performances of the 
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NAG solver. Any change of these parameters should be implemented only after a careful reading of the NAG manual 

m. 

Further relevant parameters are zin and zend, setting respectively the initial and final value of the independent 
variable z. Their present values correspond to the two temperatures of T,; = lOMeV and Tf = 1/130 MeV. Note that 
a few settings depend on these values, which then should be varied with caution. 

C. INIT 

Besides initializing the nuclear parameters, this subroutine calculates the initial values for all nuclide abundances 
and the electron chemical potential, the latter requiring the inversion of the implicit equation (|21| with a NAG routine. 

D. FCN, THERMO , RATE, EQSLIN 

The subroutine FCN is required by the NAG solver to calculate the right hand side of the differential BBN equations. 
In order to do this, the thermodynamical quantities which appear in the equations are evaluated with a call to the 
subroutine THERMO. The second step is the calculation of the reaction rates with the subroutine RATE. Then the 
linearization of the set of equations is performed, with the construction of a corresponding N nuc x N nuc matrix (N nuc 
being the number of nuclides). In this way, the unknown functions appear in a linear equation system, solved by 
Gaussian elimination in the subroutine EQSLIN (this method is very similar to the one applied in the Kawano code 
0). 

E. 0UTEV0L, OUTEND 

The subroutine OUTEVOL is called during the evolution, for printing the intermediate values of the chosen nuclide 
abundances in one of the output files. Moreover, if requested by the user, this subroutine allows to follow the 
resolution evolution, printing some physical quantities on the screen. Finally, OUTEND prints the final values of the 
nuclide abundances and electron chemical potential in the other output file along with some technical information on 
the differential evolution resolution. The final yield of the i— th nuclide is expressed as the ratio Xi/X p , i.e. number 
density normalized to hydrogen. The only exceptions are Hydrogen expressed as X p and 4 He, which is conventionally 
reported in terms of the (approximate) mass fraction Y p — AX± Ho . 

V. MAIN DIFFERENCES WITH RESPECT TO THE WAGONER-KAWANO CODE 

As we have previously emphasized, the public Kawano code [20j ] was the starting point for the development of 
PArthENoPE . So, they have a similar structure, like the subdivision in several subroutines which contain, for example, 
the interface with the user, the calculation of nuclear rates, of the thermodynamical quantities, the differential equation 
solver, and the production of the output. Here we summarize the main physical and numerical differences with respect 
to the original code: 

• The interface menu and the options are now different. In particular, the user can more easily implement changes 
in the nuclear reaction network. A card-mode input is available. The output is easily customized. 

• Several numerical routines have been replaced with more efficient algorithms, most of which using NAG routines. 

• Improved calculations for the n — p reactions are implemented via new fits, not as effective corrections added a 
posteriori. They also include effects of finite nucleon mass and non-thermal neutrino spectral distortions. The 
same holds for the case with neutrino asymmetry. 

• An improved and more accurate calculation of thermodynamical variables is implemented, both in the elec- 
tromagnetic sector and the neutrino one. In particular, for the latter case we do not simply impose entropy 
conservation, but take into account entropy transfer during the e + — e~ annihilation phase. 

• We implemented an updated nuclear reaction network, obtained using new data and reduction techniques. Also, 
new reactions and new types of reactions (different number of nuclides in the initial and final states) are now 
included. The new code has been restructured so that it is easier to implement new reactions and that only 
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FIG. 2: The relative difference in percent of nuclei abundances between PArthENoPE and the Kawano code, AXi/Xi = 2(Xf — 
Xf < )/(Xf versus baryon density Os/i 2 . The solid (dashed) curves are for positive (negative) values. 

fundamental nuclear data are needed as input. All derived quantities are calculated in the code. This is aimed 
at simplifying future updates. 

• We improved the numerical resolution of the coupled BBN equations, by using a multistep method, belonging to 
the class of Backward Differentiation Formulas (instead of the traditional Runge-Kutta solver of the Kawano code 
which is a single step-method) implemented by a NAG routine. The default values of the accuracy parameters 
have been chosen to guarantee a good compromise between the accuracy goal and a reasonable running time. In 
particular the relative accuracy reached on 4 He mass fraction is of the order of 10 -4 , thus keeping the numerical 
error below the level of theoretical uncertainties. 

In Figure[5]we show the relative variation (in percent) of nuclear abundances computed with PArthENoPE and using 
the original Kawano code vs. baryon density, for the fiducial values of all the other parameters. To give an idea 
of the differences introduced by accounting for the major updates to the physics implemented in the code, in Table 
IVIII we report, for Qgh 2 = 0.022 (i.e. consistent with the value singled out by the WMAP Collaboration [l[), the 
relative variation with respect to the Kawano code results of the theoretical prediction of abundances in three cases: 
1) improved treatment of n-p weak rates, as presently in PArthENoPE, but original nuclear network as in Kawano code 
(first column); 2) improved treatment of n-p weak rates as presently in PArthENoPE and updated nuclear network as 
in Ref. Q (second column); 3) PArthENoPE code, with improved calculation of thermodynamical variables, both in 
the electromagnetic sector and the neutrino one, and complete improvement of the nuclear network (third column). 
In summary, the PArthENoPE prediction for 4 He differs with respect to the original Kawano code by 0.7 %, mainly 
due to the improved treatment of radiative corrections, finite temperature and finite nucleon mass corrections in 
neutron-proton weak rates. Actually, this value is larger than the theoretical accuracy of PArthENoPE on Y p , which 
is of the order of 0.2%, see also @. On the other hand, our results are, for example, in very good agreement with 
[43| , where all mentioned corrections to weak rates are included as discussed in details in [35[ . Rescaling the value of 
neutron lifetime to the values r„ = 886.7 s adopted in [43| and using their fit for Y p the agreement is at the 0.1% level. 
Concerning the other nuclides, the variation of the theoretical predictions with respect to the Kawano code reaches 
larger values, of the order of 3 % for 6 Li, at the 1% level for 7 Li and 7 Be, while it is very small for 2 H, as in this case 










9 Be/H 

10 g /H 




\ "C/H 


: s l i /H 
ui/ n ...... 




■ 3 He/H . 

" 7 Li/H : >\ 

" 7 Be/H ^ \ 


lie/ ti 

. - 7 Li/H 


- 4He — V^V/^/ 




- 2 H/H \ \ / 
3 H/H f , \/ , A/ 


2 H/H 
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Nuclide 


weak rates (%) 


nucl. rates (%) 


Parthenope (%) 


2 H 


1.5 


3.2 


0.2 


3 H 


3.3 


2.3 


-0.7 


3 He 


0.1 


-3.4 


-1.1 


4 He 


0.4 


0.4 


0.7 


6 Li 


-2.7 


12.2 


2.8 


7 Li 


-1.8 


2.9 


1.3 


7 Be 


-2.0 


2.5 


1.6 


9 Be 


-0.1 


200.0 


200.0 


io B 


1.3 


-150.0 


-153.8 


U C 


200.0 


146.0 


153.8 


!2 C 


0.4 


0.5 


-43.6 



TABLE VII: The relative variation in percent of the theoretical predictions between PArthENoPE and the Kawano code, 
AX Z /X Z = 2(Xf - Xf)/(Xf + Xf), is shown for three cases (see text). The results are shown for Q, B h 2 = 0.022 and 
standard values for all other physical parameters. 

the effects of improvements on weak rates treatment and updated nuclear rate network have different sign with respect 
to the one in plasma and neutrino treatment and almost cancel out accidentally. We notice that the introduction of 
new processes (reactions 98, 99, 100 of Table H*V)) . see [22J, and the update of those already considered in the Kawano 
code result in a large difference in the BBN theoretical prediction for metallicity, although probably insufficient to 
change the chemistry of primordial clouds. Finally, the reader may want to consider similar comparisons performed 
in the literature between updated versions of BBN codes used, e.g. in [44[, [45[ and the already quoted [43|, and the 
results of the Kawano code presented in [2l| . 

VI. CONCLUSIONS 

In this paper we have described the general structure and features of PArthENoPE , a new numerical code which 
computes the theoretical abundances of nuclei produced during BBN, as function of several input cosmological param- 
eters. This code has been recently made public and can be obtained at the URL http://parthenope.na.infn.it/. The 
code evaluates the abundances of 26 nuclide in the standard BBN scenario, as well as in extended models allowing 
for extra relativistic particles or neutrino chemical potential. We checked for a few fiducial cases that the results for 
the 4 He abundance agree very closely with those of Ref. [43j , where the appropriate corrections [35| are applied to 
the original Wagoner/Kawano code |19l [20l|. 

In view of the improved data coming from astrophysical observations, accurate tools providing theoretical predictions 
on cosmological observables are required to check the overall consistency of the picture of the evolution of the Universe, 
as well as for investigating and constraining new physics beyond the present framework of fundamental interactions. 

Much effort has been put in the recent years by several groups in order to increase the level of accuracy of theoretical 
prediction on nuclide abundances, in particular by improving the estimate of the neutron to proton weak conversion 
rates and the nuclear network rates. The results of these studies, along with a to date analysis of experimental results 
on relevant nuclear reactions have been implemented in PArthENoPE , which hopefully will turn useful as an accurate 
tool for BBN-related studies. 

Acknowledgments 

In Naples, this work was supported in part by the PRIN04 "Fisica Astroparticellare e Cosmologia" and PRIN06 
"Fisica Astroparticellare: neutrini ed universo primordiale" of the Italian MiUR. P.D.S. acknowledges support by the 
US Department of Energy and by NASA grant NAG 5- 10842. 



13 



APPENDIX A: DERIVATION OF THE PARTHENOPE SET OF EQUATIONS 

We define z = m e /T, x — m e a, z — x/z = aT = T/T u , Ub = m^ 3 ns and introduce the following quantities: 



=x(z) 



, p v = a p v = — p v 



P = Pe-yB + Pv , P = Pe 7 B + Pz, , 



Starting from Eq.s © and ©, 



3if , 



n B 

TLB 

P= -3i? (p + p) , 



and separating the neutrino contribution one gets 

Pe- t B + Pv = -3 if {p e ~,B + Pe 7 s) - 4 if /5„ , 

where in (|A7[) we have used p v — ip v . From Eq. (| A6|) one gets the time derivative 



rn t . 



P=\ — ) P + 4 — — P 



and thus 



Pe-yB = - 3 if {pe~ ( B + Pe-fB ) + 4 -p e7 B 



2 



For the neutrino energy density one gets 



dpv . . dp v dp v 

p v = — — x — m e a— — = tlx— — 

ax ax ax 



if 



X ~T^ ~ ^Pv 

ax 



. rr /"m e \ 4 Tr dp v 

p v + AH Pu = — Hx^ . 

\ x ) ax 

Hence substituting (|A10|) in (IA9|) we obtain 

Z 

Pe-yB = 4 - Pe-yB - 3 if {p el B + pVyfi) - H Af(z) 



(Al) 

(A2) 
(A3) 
(A4) 

(A5) 
(A6) 

(A7) 
(A8) 

(A9) 



(A10) 



(All) 



On the other hand, the total time derivative of p e -yB can be expressed via the partial derivatives with respect to z, 
(j> e and Xi , 



d Pe-yB . , dp el B 
Pe-yB = — z + 



dz 



d(j)e 



E 



dp, 



e-yB 



8X t 



Xi 



dpe-jB dpejB 



dz 



d(f)e dz 



d(p e \ - dpe-yB dX % \ . 

J* dXi dz i Z 



Thus, equating the r.h.s of (lAlip and (|A12|) after some rearrangement reads 



dp ei B 4 



dpe-yB d(j} e 

dz z rL7 ° ' Q(h d z ' jL-j QXi dz 



Pe-fB 



E 



dpe-yB dXi 



3zH (pe-yB + Pe-yB) - H Af(z) 



Starting from (|A5|) and proceeding in the same way leads to 



dfis OflB d<fi e ®?ib dXi 

dz 84> e dz t-j QXi dz 



z = — 3 z if its 



(A12) 



(A13) 



(A14) 
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Obtaining z from (|A14|) and substituting into (|A13|) we obtain 

:PejB 



3 H 7%b 



dpe-,B 4 » | dPe-fB d<j>^ , dp c ~,B (IXj 

d<p c dz ~*~ OX; dz 



Wb i gng dffle , oris dXj 

dz ' dcj> e dz ' l^i dXi dz 



By using Eq. © 



n B 



L{z,<j)e) 



= -3zH {p eiB + b e7B ) - HAf(z) . (A15) 



(A16) 



one can express fiB and its derivatives as function of L(z,cf) e ) which is defined in (jlip . By solving Eq. (|A15[) with 
respect to d<j> e /dz one gets 



U(z) 



d(p e I LKi + [ p n B + Pe-yB + ) '<'-> 

dz z 



dp? dt ( ~ i - , AT(z) 



where 



(A17) 



-^• + « + ^--4-=t4( 3i -4) & -E^?( 4ii+ s) f - < Ai8 » 



K2 = z — 3 L — z L 

oz 



(A19) 



According to our notations, AM S and M„ stand for the i-th nuclide mass excess and the atomic mass unit, respectively, 
normalized to m e , whereas H = m e H and I\ = m e Ti. By substituting in (|A14[) the expression obtained for d<fi e /dz 
in (|A17|) we get 



oh 



ZH 



77~ V Pe 7 S + Pe-yB + -3- j ~ n B -g£- 



9p c 



dfiB ( dPi-tB _ 4 ~ , dpe-yB dXj \ _ dfie-iB ( dh B 1 dh B dXj 

0^ \ dz z reyB T Z^i dz J 90 c \^ dz + Z^i ax s dz 



namely 



9L 



3i? 



Pe-yB + Pe~fB 



AT(z) 



The equations for the abundances ([7]) then become 

dXj ■ dt r, 



xA = 

3zH L- 



^Pe 7 B + Pe 7 B + -3- 



The solution of neutrino dynamics performed in [131 Il4 , Il5| | allows to compute the quantity 



1 



(A20) 



(A21) 



(A22) 



(A23) 



Notice that N{z) would vanish for purely thermal neutrinos, and it is strictly related to the small entropy transfer to 
neutrinos during the e + — e~ annihilation stage. In the code, N(z) is calculated by using the following fit, which is 
accurate to better than 0.1% in the relevant range: 



ex P fe/=i n i 2 ') z < 4 > 
z > 4 , 



(A24) 



10.21703221236002 n x - 
-340.3323864212157 n 3 = 
-2045.577491331372 n 5 = 
2266.1521815470196 n 7 - 
586.0618273295763 n 9 = 
-35.715878215468045 n n 
-0.3713438862054167 n 13 



= 61.24438067531452 
= 1057.2707914654834 
= 2605.9087171012848 
= 1374.2623075963388 
= 174.87532902234145 
= 4.7538967685808755 
= 0.012908416591272199 
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APPENDIX B: COMMON VARIABLES USED IN PARTHENOPE 



VARIABLE 


DESCRIPTION 


COMMON 


A A /attvttt \ 

AA(NNUC) 


Nuclide atomic numbers, A4 


A TV TT TK IT 

ANUM 


1 AC I OR(NREC) 


Multiplicative lactor lor the rate of reaction l-tn 


CHRA1E 


HCHRAT(NREC) 


Type of changes adopted for reaction i-th 




NCHRAT 


Number of reactions to be changed 




WCHRAT(NREC) 


Reactions to be changed 




ALF 


Fine structure constant, a 


CONSTANTS 


COEF(4) 


Unity conversion factors 




GN 


Newton constant, Gat, in MeV~ 2 




ME 


Electron mass, m e 




MU 


Atomic mass unit, M u 




PI 


TT 




IFCN 


Counter 


COUNTS 


IFCN1 


Counter 




IS AVE 1 


Counter 




ISAVE2 


Counter 




TO T^T^Tl 

1S1EP 


Counter 




DZ 


Stepsize ol the independent variable m the resolution ol the 
nucleosynthesis equations 


TATT'T A r / 

DELIA/ 


DZO 


Initial value for DZ 




TAA J ( C\. AT ATT T/~*\ 

DM(U:JNJN UOj 


Mass excesses in MeV, AM S 


UMAbb 


DMfi(NNUC) 


A 1 • '1 Air/ 

Adimensional mass excesses, i\Mi/m e 


TA1V T A PCTT 

DMASSH 


PHI 


Adimensional electron chemical potential, <p e 


ECHPO I 


A /^i / TV TTV TT T /^\ A \ 

AG(NNUC,4) 


"XT 1 1 * i * f* 1 • fT* • i 

INuclear partition lunction coefficients 


GPART 


YY0(NNUC+1) 


t*j*i i rij i *i i 1*11 l i*i 

Initial values ol electron chemical potential, <p e , and nuclide 
abundances, X; 


"TTV T A 1 1 1 TXT 

INABUN 


CMODE 


Tii r j ~] i * r j i • l 

l 1 lag lor the choice ol the running mode 


INPCARD 


FOLLOW 


Option lor lollowmg the evolution on the screen (card mode] 




OVERW 


Option for overwriting the output files (card mode) 




INC 


Maximum value ol the nag lor the convergence ol the matrix inversion 


INVFLAGS 


MBAD 


Error flag for the matrix inversion 




LHO 


Initial value lor the adimensional electron/ positron asymmetry, L 


INVPHI 


zo 


Initial value lor the evolution variable z = m e /l (=ZIN) 




AMAT(NNUC,NNUC) 


TVTj* • 1 1' i 1 1 • • j • p i 1 1 j • 1 j \ ■ / i 7 \ 

Matrix involved m the linearization ol the relation between Xi(z + az) 
and Xi(z) 


LINCOEF 


BVEC(NNUC) 


Vector involved in the linearization of the relation between X, L {z — dz) 
and Xi(z) (contains A; m reverse order) 




~\7"V / AT ATT T/^A 

YX(NNUC) 


Xi(z) m reverse order 




YMIN 


Numerical zero of nuclide abundances 


MINABUN 


DNNU 


Number of extra effective neutrinos, AJV e // 


MODPAR 


DNNUXI 


Contribution to the number of extra effective neutrinos from a non 
zero neutrino chemical potential, AA^/ / of Eq. (|25|) 




ETAF 


■ »t 1 pjii j ij l "i i' 

Imal value ol the baryon to photon density ratio, r/j 




IXIE 


A • j * *j C* * i l 1j j* l • l j j • 1 

A positive integer hxing the electron neutrino chemical potential 




RHOLMBD 


bncrgy density corresponding to a cosmological constant, p\ 




TAU 


Value of neutron lifetime in seconds, r n 




XIE 


T^leftron neutrino ehpTnie/d notenti/d P ( — XTT^fM 




XIEO(NXIE) 


Electron neutrino chemical potential, ^ 


NCHPOT 


INUC 


Number of nuclides in the selected network 


NETWRK 


IREC 


Number of reactions among nuclides in the selected network 




IXT(30) 


Code of the nuclides whose evolution has to be followed 






(ixt(30)=control integer) 





17 



VARIABLE 


DESCRIPTION 


COMMON 


BYY(NNUC+1) 


Text strings for the output 


NSYMB 


MN 


Neutron mass, M n 


NUCMASS 


MP 


Proton mass, M p 




NAMEFILE1 


Name of the output file for the final values of the nuclide abundances 


OUTFILES 


NAMEFILE2 


Name of the output file for the evolution of the nuclides whose 
evolution has to be followed 




NBH 


Adimensional baryon number density, n^/m^ 


OUTVAR 


THETAH 


Adimensional Hubble function times 3, 3H 




rrvTT 

1XH 


INeutrmo to pnoton temperature ratio, ijf/J 




U/jr 


Previous iteration value of the step-size of the independent variable 
in the resolution of the nucleosynthesis equations 


1 i v I ^ V V Jr\ 1 . 


bUMMYr 


Value of the linear combination 1 AM; + ^ 1 1 i 




SUMZYP 


Value of the linear combination Zi Yi 




ZP 


Previous iteration value of the evolution variable z = m e /T 




CFLAG 


Flag for the input variable type in the card reading (card mode) 


READINP 


IKEY 


Progressive argument key number in the card reading (card mode) 




ISTART 


Starting point of the line in the card reading (card mode) 




DNCHRAT 


Number of reactions to be added to the changed ones 




LINE 


Line input from the card file (card mode) 




IFORM(NREC) 


Reaction type (1-12) 


RECPAR 


NG(NREC) 


Number of incoming nuclides of type TG 




NH(NREC) 


Number of incoming nuclides of type TH 




NI(NREC) 


Number of incoming nuclides of type TI 




NJ(NREC) 


Number of incoming nuclides of type T J 




NK(NREC) 


Number of incoming nuclides of type TK 




NL(NREC) 


Number of incoming nuclides of type TL 




Q9(NREC) 


Energy released in reaction (in unit of 10 9 K) 




REV(NREC) 


Reverse reaction coefficient 




TG(NREC) 


Incoming nuclide type 




TH(NREC) 


Outgoing nuclide type 




TI(NREC) 


Incoming nuclide type 




TJ(NREC) 


Incoming nuclide type 




TK(NREC) 


Outgoing nuclide type 




TL(NREC) 


Outgoing nuclide type 




RATEPAR(NREC,13) 


Reaction parameter values (=IFORM+TI+...+NI+...) 


RECPARO 


RSTRING(NREC) 


Reaction text strings 


RSTRINGS 


LH 


Function L 


THERMQ 


LHPHI 


Derivative of L with respect to (f> e 




LHZ 


Derivative of L with respect to z 




NAUX 


Neutrino auxiliary function, N{z) 




PBH 


Adimensional baryon pressure, ps 




PEH 


Adimensional electron pressure, p e 




PGH 


Adimensional gamma pressure, p 7 




RHOBH 


Adimensional baryon energy density, ps 




RHOEH 


Adimensional electron energy density, p e 




RHOEHPHI 


Derivative of p e with respect to cf) e 




RHOEHZ 


Derivative of p e with respect to z 




RHOGH 


Adimensional gamma energy density, p 7 




RHOGHZ 


Derivative of /5 7 with respect to z 




RHOH 


Adimensional total energy density, p 




GNUC(NNUC) 


Nuclide spin degrees of freedom 


SPINDF 
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VARIABLE 


DESCRIPTION 


COMMON 


A(13) 
B(10) 

DA(12,NXIE) 

DB(12,NXIE) 

QNP 

QNPl 

QPN 

QPNl 


Forward weak reaction best-fit parameters (non degenerate case) 
Reverse weak reaction best-fit parameters (non degenerate case) 
Forward weak reaction best-fit parameters (degenerate case) 
Reverse weak reaction best-fit parameters (degenerate case) 
Forward weak reaction best-fit exponent parameter 
Forward weak reaction best-fit exponent parameter 
Reverse weak reaction best-fit exponent parameter 
Reverse weak reaction best-fit exponent parameter 


WEAKRATE 


ZZ(NNUC) 


Nuclide atomic charges, Zi 


ZNUM 
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